#ifndef FCTSPACE_
#define FCTSPACE_

#include <iostream>
#include <fstream>
#include <sstream>
#include <string>
#include <vector>
#include <map>
#include <cmath>

#include "DGMatrix.h"
#include "DGIntMatrix.h"
#include "ANN/ANN.h"
#include "Parser/Parser.h"

using namespace std;

enum FctType {LagrangeTri=2,LagrangeQuad=3};
//------------------------------------------------------------------------------
class FunctionSpace{

private:
 
public:
  
  FunctionSpace(){};
  ~FunctionSpace(){};
  
  virtual void setQuadratureRules ()=0;
  virtual void setRefNodes        ()=0;
  virtual void setMatrices        ()=0;
  
  virtual int                 getOrder            ()=0;  
  virtual int                 getNbENP            ()=0;
  virtual int                 getNbEQP            ()=0;
  virtual int                 getNbFNP            ()=0;
  virtual int                 getNbFQP            ()=0;  
  
  virtual double              getERefNodesPos     (int,int)=0;
  virtual const DGMatrix*     getERefNodesPos     ()=0;
  virtual int                 getFRefNodesIndex   (int,int)=0;
  virtual const DGIntMatrix*  getFRefNodesIndex   ()=0;
  virtual double              getERefNodesCoef    (int,int)=0;
  virtual const DGMatrix*     getERefNodesCoef    ()=0;
  virtual double              getFRefNodesCoef    (int,int)=0;
  virtual const DGMatrix*     getFRefNodesCoef    ()=0;
  virtual double              getEpp              (int,int)=0;
  virtual const DGMatrix*     getEpp              ()=0;
  
  virtual double              getEQP              (int,int)=0;
  virtual const DGMatrix*     getEQP              ()=0;
  virtual double              getFQP              (int)=0;
  virtual const DGVector*     getFQP              ()=0;
  virtual double              getEW               (int)=0; 
  virtual const DGVector*     getEW               ()=0; 
  virtual double              getFW               (int)=0;
  virtual const DGVector*     getFW               ()=0; 
  
  virtual const DGMatrix*     getPhiEQP           ()=0;
  virtual const DGMatrix*     getPhiFQP           ()=0;
  virtual const DGMatrix*     getdPhiUVPhiEQP     (int,int)=0;
  virtual const DGMatrix*     getPhiPhiFQP        (int)=0;
  
  virtual const DGMatrix*     getFqduv            (int)=0;
  virtual const vector<DGMatrix>* getFqduv        ()=0;
  virtual const DGMatrix*     getFquv             ()=0;
  
  virtual const DGMatrix*     getMass             ()=0;
  virtual const DGMatrix*     getInvMass          ()=0;
  virtual const DGMatrix*     getLumpedMass       ()=0;
  virtual const DGMatrix*     getInvLumpedMass    ()=0;
  
  virtual double              getPhiF             (int,double)=0;
  virtual double              getPhiE             (int,double,double)=0;
  virtual double              getdPhiE            (int,int,double,double)=0;

};
//------------------------------------------------------------------------------
class LagrangeFctSpaceForTriangle:public FunctionSpace{

public:
  
  string FSname;
  
  LagrangeFctSpaceForTriangle(int);
  ~LagrangeFctSpaceForTriangle(){};
  
private: 
  
  int order;
  int Qorder;
  int nbEQP;
  int nbFQP;
  int nbENP;
  int nbFNP;
  
  DGMatrix EQP;                                                                 // (nbEQP,dim)
  DGVector EW;                                                                  // (nbEQP)
  DGVector FQP;                                                                 // (nbFQP)
  DGVector FW;                                                                  // (nbFQP)  
  
  DGMatrix Mass;                                                                // (nbENP,nbENP)
  DGMatrix InvMass;
  DGMatrix LumpedMass;
  DGMatrix InvLumpedMass;
  
  DGMatrix    ERefNodesPos;                                                     // (nbENP,dim)
  DGMatrix    ERefNodesCoef;                                                    // (nbENP,nbENP)
  DGMatrix    FRefNodesCoef;                                                    // (nbFNP,nbFNP)
  DGIntMatrix FRefNodesIndex;                                                   // (nbFaces,nbFNP)
  DGMatrix    Epp;                                                              // (nbENP,dim)
  
  DGMatrix phiEQP;                                                              // (nbEQP,nbENP)
  DGMatrix phiFQP;                                                              // (nbFQP,nbFNP)
  
  vector<vector<DGMatrix> > dphiUVphiEQP;                                       // [dim][i](nbENP,nbEQP)
  vector<DGMatrix>          phiphiFQP;                                          //      [i](nbFNP,nbFQP)

  vector<DGMatrix> dphiUVEQP;                                                   // [dim](nbEQP,nbENP)
  vector<DGMatrix> Fqduv;                                                       // [dim](nbENP,nbEQP)
  DGMatrix         Fquv;                                                        // (nbFNP,nbFQP)
  
public:
  
  double monomeof   (double,double,int);
  double dmonomeof  (double,double,int,int);
  
  virtual void setQuadratureRules ();
  virtual void setRefNodes        ();
  virtual void setMatrices        ();
  
  virtual int                 getOrder            ()                    {return order;};  
  virtual int                 getNbENP            ()                    {return nbENP;}
  virtual int                 getNbEQP            ()                    {return nbEQP;}
  virtual int                 getNbFNP            ()                    {return nbFNP;}
  virtual int                 getNbFQP            ()                    {return nbFQP;}  
  
  virtual double              getERefNodesPos     (int nI,int dim)      {return ERefNodesPos(nI,dim);}
  virtual const DGMatrix*     getERefNodesPos     ()                    {return &ERefNodesPos;}
  virtual int                 getFRefNodesIndex   (int fI,int nI)       {return FRefNodesIndex(fI,nI);}
  virtual const DGIntMatrix*  getFRefNodesIndex   ()                    {return &FRefNodesIndex;}
  virtual double              getERefNodesCoef    (int nI_i,int nI_j)   {return ERefNodesCoef(nI_i,nI_j);}
  virtual const DGMatrix*     getERefNodesCoef    ()                    {return &ERefNodesCoef;}
  virtual double              getFRefNodesCoef    (int nI_i,int nI_j)   {return FRefNodesCoef(nI_i,nI_j);}
  virtual const DGMatrix*     getFRefNodesCoef    ()                    {return &FRefNodesCoef;}
  virtual double              getEpp              (int nI_i,int dim)    {return Epp(nI_i,dim);}
  virtual const DGMatrix*     getEpp              ()                    {return &Epp;}
  
  virtual double              getEQP              (int qI,int dim)      {return EQP(qI,dim);}
  virtual const DGMatrix*     getEQP              ()                    {return &EQP;}
  virtual double              getFQP              (int qI)              {return FQP(qI);}
  virtual const DGVector*     getFQP              ()                    {return &FQP;}
  virtual double              getEW               (int wI)              {return EW(wI);} 
  virtual const DGVector*     getEW               ()                    {return &EW;} 
  virtual double              getFW               (int wI)              {return FW(wI);} 
  virtual const DGVector*     getFW               ()                    {return &FW;} 

  virtual const DGMatrix*     getPhiEQP           ()                    {return &phiEQP;}
  virtual const DGMatrix*     getPhiFQP           ()                    {return &phiFQP;}
  virtual const DGMatrix*     getdPhiUVPhiEQP     (int dim,int i)       {return &dphiUVphiEQP[dim][i];}
  virtual const DGMatrix*     getPhiPhiFQP        (int i)               {return &phiphiFQP[i];}

  virtual const DGMatrix*     getFqduv            (int dim)             {return &Fqduv[dim];}
  virtual const vector<DGMatrix>* getFqduv        ()                    {return &Fqduv;}
  virtual const DGMatrix*     getFquv             ()                    {return &Fquv;}

  virtual const DGMatrix*     getMass             ()                    {return &Mass;};
  virtual const DGMatrix*     getInvMass          ()                    {return &InvMass;};
  virtual const DGMatrix*     getLumpedMass       ()                    {return &LumpedMass;};
  virtual const DGMatrix*     getInvLumpedMass    ()                    {return &InvLumpedMass;};

  virtual double              getPhiF             (int,double);
  virtual double              getPhiE             (int,double,double);
  virtual double              getdPhiE            (int,int,double,double);  
  
  
};
//------------------------------------------------------------------------------
class FunctionSpaceFactory{
public:
  static FunctionSpace* CreateNewFunctionSpace(FctType,int);
};

#endif
